2024-04-03
#install.packages("~/Desktop/tutorial/space_0.1-1.1.tar.gz", repos = NULL, type = "source")
#install_github(repo = "topherconley/spacemap")
suppressPackageStartupMessages(library(igraph))
suppressPackageStartupMessages(library(ggplot2))
suppressPackageStartupMessages(library(ggpubr))
suppressPackageStartupMessages(library(space))
suppressPackageStartupMessages(library(glasso))
#suppressPackageStartupMessages(library(spacemap))Estimates a sparse inverse covariance matrix using a lasso (L1) penalty
glasso(s, rho, nobs=NULL, zero=NULL, thr=1.0e-4, maxit=1e4, approx=FALSE,
penalize.diagonal=TRUE, start=c("cold","warm"), w.init=NULL,wi.init=NULL, trace=FALSE)
s : Covariance matrix:p by p matrix (symmetric)
rho : (Non-negative) regularization parameter for lasso. rho=0 means no regularization.
data(spaceSimu)
Y.m<- spaceSimu$Y.data
p<- dim(Y.m)[2] # no. of nodes
n=nrow(Y.m)
s<- var(Y.m)
a<-glasso(s, rho=0.25) ####### change the rho parameter to control the number of edges.
a$w[1:10,1:10]## [,1] [,2] [,3] [,4] [,5]
## [1,] 1.250000e+00 -7.133126e-03 6.382393e-03 9.811987e-11 2.015679e-04
## [2,] -7.133126e-03 1.250000e+00 -4.257136e-05 -1.726507e-08 -1.344484e-06
## [3,] 6.382393e-03 -4.257136e-05 1.250000e+00 -6.239176e-13 3.947731e-02
## [4,] 9.811987e-11 -1.726507e-08 -6.239176e-13 1.250000e+00 -4.957146e-14
## [5,] 2.015679e-04 -1.344484e-06 3.947731e-02 -4.957146e-14 1.250000e+00
## [6,] 9.372103e-03 -5.348190e-05 4.785325e-05 6.607423e-13 1.511298e-06
## [7,] 1.042896e-03 -6.956248e-06 2.042524e-01 -2.564615e-13 2.415964e-01
## [8,] 1.995857e-05 -1.138977e-07 1.019207e-07 3.229361e-14 3.219189e-09
## [9,] -6.325406e-04 1.108456e-01 -3.775078e-06 -1.531005e-09 -1.192241e-07
## [10,] 3.688676e-03 -2.191352e-05 1.926529e-05 2.938483e-13 6.084340e-07
## [,6] [,7] [,8] [,9] [,10]
## [1,] 9.372103e-03 1.042896e-03 1.995857e-05 -6.325406e-04 3.688676e-03
## [2,] -5.348190e-05 -6.956248e-06 -1.138977e-07 1.108456e-01 -2.191352e-05
## [3,] 4.785325e-05 2.042524e-01 1.019207e-07 -3.775078e-06 1.926529e-05
## [4,] 6.607423e-13 -2.564615e-13 3.229361e-14 -1.531005e-09 2.938483e-13
## [5,] 1.511298e-06 2.415964e-01 3.219189e-09 -1.192241e-07 6.084340e-07
## [6,] 1.250000e+00 7.819330e-06 1.496433e-07 -4.742587e-06 2.765652e-05
## [7,] 7.819330e-06 1.250000e+00 1.665582e-08 -6.168556e-07 3.147987e-06
## [8,] 1.496433e-07 1.665582e-08 1.250000e+00 -1.010005e-08 5.889660e-08
## [9,] -4.742587e-06 -6.168556e-07 -1.010005e-08 1.250000e+00 -1.943214e-06
## [10,] 2.765652e-05 3.147987e-06 5.889660e-08 -1.943214e-06 1.250000e+00
## [,1] [,2] [,3] [,4] [,5] [,6]
## [1,] 1.007181838 0.004684408 0.0000000 0.000000 0.0000000 -0.005998483
## [2,] 0.004684408 0.849759746 0.0000000 0.000000 0.0000000 0.000000000
## [3,] 0.000000000 0.000000000 0.9017998 0.000000 0.0000000 0.000000000
## [4,] 0.000000000 0.000000000 0.0000000 0.800471 0.0000000 0.000000000
## [5,] 0.000000000 0.000000000 0.0000000 0.000000 0.8310445 0.000000000
## [6,] -0.005998483 0.000000000 0.0000000 0.000000 0.0000000 0.861787423
## [7,] 0.000000000 0.000000000 -0.1286478 0.000000 -0.1606219 0.000000000
## [8,] 0.000000000 0.000000000 0.0000000 0.000000 0.0000000 0.000000000
## [9,] 0.000000000 -0.071502661 0.0000000 0.000000 0.0000000 0.000000000
## [10,] 0.000000000 0.000000000 0.0000000 0.000000 0.0000000 0.000000000
## [,7] [,8] [,9] [,10]
## [1,] 0.0000000 0.0000000 0.00000000 0.0000000
## [2,] 0.0000000 0.0000000 -0.07150266 0.0000000
## [3,] -0.1286478 0.0000000 0.00000000 0.0000000
## [4,] 0.0000000 0.0000000 0.00000000 0.0000000
## [5,] -0.1606219 0.0000000 0.00000000 0.0000000
## [6,] 0.0000000 0.0000000 0.00000000 0.0000000
## [7,] 0.8572223 0.0000000 0.00000000 0.0000000
## [8,] 0.0000000 0.8304466 0.00000000 0.0000000
## [9,] 0.0000000 0.0000000 0.80634068 0.0000000
## [10,] 0.0000000 0.0000000 0.00000000 0.8186895
## [1] 781
temp=graph.adjacency(adjmatrix=fit.adj, mode="undirected")
#plot(temp, vertex.size=3, vertex.frame.color="white",layout=layout.fruchterman.reingold, vertex.label=NA, edge.color=grey(0.5))Degree of every node = number of edges this node has with other nodes.
Hub nodes: Nodes connected with large number of other nodes.
temp.degree=apply(fit.adj, 2, sum)
tab<- cbind.data.frame("Nodes"=colnames(Y.m), "degree"=temp.degree)
s<- sort(tab$degree, ind=T, decreasing = T)
tab.sort<- tab[s$ix,]
tab.sort$Nodes<- factor(tab.sort$Nodes, levels=unique(tab.sort$Nodes))
dat<- tab.sort[1:30,]
#ggplot(dat, aes(y=degree, x=Nodes)) +geom_bar(position="stack", stat="identity") + theme_bw() + theme(axis.text.x = element_text(colour="black",size=15,angle=45,hjust=1,vjust=1,face="plain"),axis.title.x = element_text(colour="black",size=25,angle=0,hjust=.5,vjust=.5,face="plain"), axis.text.y = element_text(colour="black",size=15,angle=0,hjust=1,vjust=0,face="plain"), axis.title.y = element_text(colour="black",size=25,angle=90,hjust=.5,vjust=.5,face="plain"),legend.title = element_text(colour="black",size=15), legend.text=element_text(colour="black",size=15), legend.direction = "horizontal", legend.position = "top") + theme(plot.margin = unit(c(1,1,1,1), "cm"))space.neighbor: A function to estimate partial correlations using the neighborhood selection approach
space.neighbor(Y.m, lam1, lam2=0, iter)
lam1: numeric value. This is the l_1 norm penalty parameeter.
lam2: numeric value. If not specified, lasso regression is used in the neighborhood selection.
data(spaceSimu)
Y.m<- spaceSimu$Y.data
p<- dim(Y.m)[2] # no. of nodes
n=nrow(Y.m)
alpha=1
l1=1/sqrt(n)*qnorm(1-alpha/(2*p^2))
iter=3
result1=space.neighbor(Y.m, lam1=l1*0.7, lam2=0)
result1$ParCor[1:10,1:10]## [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
## [1,] 1 0.0000000 0.00000 0 0.0000000 0 0.0000000 0 0.0000000 0
## [2,] 0 1.0000000 0.00000 0 0.0000000 0 0.0000000 0 0.1049529 0
## [3,] 0 0.0000000 1.00000 0 0.0000000 0 0.1543100 0 0.0000000 0
## [4,] 0 0.0000000 0.00000 1 0.0000000 0 0.0000000 0 0.0000000 0
## [5,] 0 0.0000000 0.00000 0 1.0000000 0 0.2620895 0 0.0000000 0
## [6,] 0 0.0000000 0.00000 0 0.0000000 1 0.0000000 0 0.0000000 0
## [7,] 0 0.0000000 0.15431 0 0.2620895 0 1.0000000 0 0.0000000 0
## [8,] 0 0.0000000 0.00000 0 0.0000000 0 0.0000000 1 0.0000000 0
## [9,] 0 0.1049529 0.00000 0 0.0000000 0 0.0000000 0 1.0000000 0
## [10,] 0 0.0000000 0.00000 0 0.0000000 0 0.0000000 0 0.0000000 1
## [1] 729
rownames(fit.adj) <- colnames(fit.adj) <- colnames(Y.m)
temp=graph.adjacency(adjmatrix=fit.adj, mode="undirected")
# plot(temp, vertex.size=3, vertex.frame.color="white",layout=layout.fruchterman.reingold, vertex.label=NA, edge.color=grey(0.5))Degree of every node = number of the parent nodes + number of child nodes
temp.degree=apply(fit.adj, 2, sum)
tab<- cbind.data.frame("Nodes"=colnames(Y.m), "degree"=temp.degree)
s<- sort(tab$degree, ind=T, decreasing = T)
tab.sort<- tab[s$ix,]
tab.sort$Nodes<- factor(tab.sort$Nodes, levels=unique(tab.sort$Nodes))
dat<- tab.sort[1:30,]
#ggplot(dat, aes(y=degree, x=Nodes)) +
# geom_bar(position="stack", stat="identity") + theme_bw() +
# theme(axis.text.x = element_text(colour="black",size=15,angle=45,hjust=1,vjust=1,face="plain"),axis.title.x = element_text(colour="black",size=25,angle=0,hjust=.5,vjust=.5,face="plain"), axis.text.y = element_text(colour="black",size=15,angle=0,hjust=1,vjust=0,face="plain"), axis.title.y = element_text(colour="black",size=25,angle=90,hjust=.5,vjust=.5,face="plain"),legend.title = element_text(colour="black",size=15), legend.text=element_text(colour="black",size=15), legend.direction = "horizontal", legend.position = "top") + theme(plot.margin = unit(c(1,1,1,1), "cm"))Modules are closely connected markers with similar functional roles/activities
obj<- graph_from_adjacency_matrix(fit.adj, mode = "undirected", weighted = TRUE)
formodule<-cluster_edge_betweenness(obj, weights = NULL, directed = FALSE, edge.betweenness = TRUE, merges = TRUE, bridges = TRUE, modularity = TRUE, membership = TRUE)
module<-communities(formodule)
size<- unlist(lapply(module,length))
mod.s5<- module[size > 5]
mod.s5[1:10]## $`1`
## [1] "g 1" "g 10" "g 11" "g 12" "g 13" "g 37" "g 46" "g 54" "g 66"
## [10] "g 99" "g 441" "g 459" "g 474" "g 480" "g 499"
##
## $`2`
## [1] "g 2" "g 19" "g 30" "g 31" "g 58" "g 59" "g 61" "g 79" "g 82" "g 83"
## [11] "g 84" "g 92"
##
## $`3`
## [1] "g 3" "g 5" "g 7" "g 14" "g 21" "g 36" "g 42" "g 49" "g 50"
## [10] "g 52" "g 63" "g 75" "g 77" "g 85" "g 90" "g 94" "g 353"
##
## $`4`
## [1] "g 4" "g 6" "g 18" "g 29" "g 38" "g 67" "g 72" "g 76" "g 95"
## [10] "g 98" "g 401" "g 408" "g 421" "g 427" "g 445"
##
## $`5`
## [1] "g 8" "g 32" "g 41" "g 43" "g 47" "g 87" "g 88" "g 411" "g 418"
## [10] "g 456" "g 468"
##
## $`6`
## [1] "g 9" "g 48" "g 93" "g 305" "g 310" "g 320" "g 332" "g 350" "g 358"
## [10] "g 361" "g 363" "g 364" "g 379" "g 389"
##
## $`7`
## [1] "g 15" "g 20" "g 25" "g 35" "g 39" "g 53" "g 60" "g 62" "g 89"
## [10] "g 409" "g 417" "g 424" "g 464" "g 484" "g 491"
##
## $`9`
## [1] "g 17" "g 45" "g 81" "g 110" "g 120" "g 143" "g 146" "g 174" "g 175"
## [10] "g 178" "g 185" "g 191" "g 198" "g 431" "g 438" "g 458" "g 471" "g 472"
## [19] "g 475" "g 482" "g 498"
##
## $`11`
## [1] "g 23" "g 28" "g 109" "g 121" "g 124" "g 135" "g 139" "g 144" "g 145"
## [10] "g 152" "g 155" "g 156" "g 164" "g 167" "g 186" "g 189" "g 194"
##
## $`15`
## [1] "g 44" "g 235" "g 307" "g 311" "g 317" "g 328" "g 335" "g 337" "g 351"
## [10] "g 354" "g 355" "g 367" "g 378" "g 384" "g 391" "g 393"
diag(fit.adj)<- F
mod<- fit.adj[which(rownames(fit.adj) %in% module[[1]]), which(colnames(fit.adj) %in% module[[1]])]
sum(mod)## [1] 30
## [1] 15 15
p.m=dim(mod)[1]
obj<- graph_from_adjacency_matrix(mod, mode = "undirected", weighted = TRUE)
# pdf("mod.spacenb.pdf",width = 40, height = 20, useDingbats=T)
# plot.igraph(obj,vertex.size=10, edge.length=5, edge.color="black", vertex.label.cex=2.5,vertex.color=V(obj)$color,layout=layout.fruchterman.reingold)
# dev.off()space.joint: A function to estimate partial correlations minimizing the joint loss function.
space.joint(Y.m, lam1, lam2=0, iter)
data(spaceSimu)
Y.m<- spaceSimu$Y.data
p<- dim(Y.m)[2] # no. of nodes
n=nrow(Y.m)
alpha=1
l1=1/sqrt(n)*qnorm(1-alpha/(2*p^2))
iter=3
result2=space.joint(spaceSimu$Y.data, lam1=l1*n*1.56, lam2=0, iter=iter)## [1] "iter=1"
## [1] "iter=2"
## [1] "iter=3"
## [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9]
## [1,] 1 0.00000000 0.0000000 0 0.0000000 0 0.0000000 0 0.00000000
## [2,] 0 1.00000000 0.0000000 0 0.0000000 0 0.0000000 0 0.09927696
## [3,] 0 0.00000000 1.0000000 0 0.0000000 0 0.1270736 0 0.00000000
## [4,] 0 0.00000000 0.0000000 1 0.0000000 0 0.0000000 0 0.00000000
## [5,] 0 0.00000000 0.0000000 0 1.0000000 0 0.2424002 0 0.00000000
## [6,] 0 0.00000000 0.0000000 0 0.0000000 1 0.0000000 0 0.00000000
## [7,] 0 0.00000000 0.1270736 0 0.2424002 0 1.0000000 0 0.00000000
## [8,] 0 0.00000000 0.0000000 0 0.0000000 0 0.0000000 1 0.00000000
## [9,] 0 0.09927696 0.0000000 0 0.0000000 0 0.0000000 0 1.00000000
## [10,] 0 0.00000000 0.0000000 0 0.0000000 0 0.0000000 0 0.00000000
## [,10]
## [1,] 0
## [2,] 0
## [3,] 0
## [4,] 0
## [5,] 0
## [6,] 0
## [7,] 0
## [8,] 0
## [9,] 0
## [10,] 1
## [1] 731
rownames(fit.adj) <- colnames(fit.adj) <- colnames(Y.m)
temp=graph.adjacency(adjmatrix=fit.adj, mode="undirected")
# plot(temp, vertex.size=3,vertex.frame.color="white",layout=layout.fruchterman.reingold, vertex.label=NA, edge.color=grey(0.5))Degree of every node = number of the parent nodes + number of child nodes
temp.degree=apply(fit.adj, 2, sum)
tab<- cbind.data.frame("Nodes"=colnames(Y.m), "degree"=temp.degree)
s<- sort(tab$degree, ind=T, decreasing = T)
tab.sort<- tab[s$ix,]
tab.sort$Nodes<- factor(tab.sort$Nodes, levels=unique(tab.sort$Nodes))
dat<- tab.sort[1:30,]
# ggplot(dat, aes(y=degree, x=Nodes)) +
# geom_bar(position="stack", stat="identity") + theme_bw() +
# theme(axis.text.x = element_text(colour="black",size=15,angle=45,hjust=1,vjust=1,face="plain"),axis.title.x = element_text(colour="black",size=25,angle=0,hjust=.5,vjust=.5,face="plain"), axis.text.y = element_text(colour="black",size=15,angle=0,hjust=1,vjust=0,face="plain"), axis.title.y = element_text(colour="black",size=25,angle=90,hjust=.5,vjust=.5,face="plain"),legend.title = element_text(colour="black",size=15), legend.text=element_text(colour="black",size=15), legend.direction = "horizontal", legend.position = "top") + theme(plot.margin = unit(c(1,1,1,1), "cm"))Modules are closely connected markers with similar functional roles/activities
obj<- graph_from_adjacency_matrix(fit.adj, mode = "undirected", weighted = TRUE)
formodule<-cluster_edge_betweenness(obj, weights = NULL, directed = FALSE, edge.betweenness = TRUE, merges = TRUE, bridges = TRUE, modularity = TRUE, membership = TRUE)
module<-communities(formodule)
size<- unlist(lapply(module,length))
mod.s5<- module[size > 5]
mod.s5[1:10]## $`1`
## [1] "g 1" "g 8" "g 10" "g 11" "g 12" "g 13" "g 16" "g 32" "g 33" "g 34"
## [11] "g 37" "g 41" "g 43" "g 46" "g 47" "g 54" "g 66" "g 80" "g 87" "g 88"
## [21] "g 99"
##
## $`2`
## [1] "g 2" "g 9" "g 17" "g 19" "g 30" "g 31" "g 45" "g 48" "g 56"
## [10] "g 58" "g 59" "g 61" "g 72" "g 79" "g 81" "g 82" "g 83" "g 84"
## [19] "g 92" "g 93" "g 100"
##
## $`3`
## [1] "g 3" "g 5" "g 7" "g 21" "g 36" "g 49" "g 50" "g 52" "g 63" "g 90"
## [11] "g 94"
##
## $`4`
## [1] "g 4" "g 29" "g 67" "g 75" "g 76" "g 85" "g 95" "g 98" "g 349"
## [10] "g 353"
##
## $`5`
## [1] "g 6" "g 18" "g 38" "g 401" "g 408" "g 421" "g 427" "g 445" "g 478"
##
## $`7`
## [1] "g 15" "g 20" "g 25" "g 39" "g 53" "g 62"
##
## $`9`
## [1] "g 23" "g 28" "g 121" "g 124" "g 139" "g 144" "g 156" "g 164" "g 189"
##
## $`14`
## [1] "g 44" "g 60" "g 307" "g 311" "g 324" "g 328" "g 355" "g 378" "g 384"
## [10] "g 391" "g 393"
##
## $`18`
## [1] "g 101" "g 105" "g 106" "g 109" "g 116" "g 122" "g 125" "g 128" "g 130"
## [10] "g 132" "g 135" "g 148" "g 152" "g 155" "g 158" "g 160" "g 163" "g 176"
## [19] "g 180" "g 184" "g 195" "g 200"
##
## $`19`
## [1] "g 102" "g 108" "g 113" "g 119" "g 127" "g 131" "g 134" "g 149" "g 150"
## [10] "g 154" "g 157" "g 168" "g 182" "g 192" "g 196"
# load("simdata.Rdata")
# p<- dim(sim1$Y)[2] # no. of nodes
# colnames(sim1$Y)<- paste0("Y_",1:p)
# net <- spacemap(Y = sim1$Y, X = sim1$X, lam1 = 70, lam2 = 28.8, lam3 = 12.38)
# net$ParCor[1:5,1:5]
# net$sig.fit[1:10]
#
# adjnet <- adjacency(net)
# dim(adjnet$yy)
# dim(adjnet$xy)
#
# adjnet$yy[1:5,1:5]
# adjnet$xy[1:5,1:5]
# adj.yy<- adjnet$yy
# adj.xy<- adjnet$xy
#
# colnames(adj.yy)<- colnames(sim1$Y)
#
# obj<- graph_from_adjacency_matrix(adj.yy, mode = "undirected", weighted = TRUE)
# pdf("netyy.pdf",width = 30, height = 30, useDingbats=T)
# plot.igraph(obj,vertex.size=10, edge.length=5, edge.color="black", vertex.label.cex=4,vertex.color=V(obj)$color,layout=layout.fruchterman.reingold)
# dev.off()# adj<- adj2igraph(adj.yy, adj.xy, yinfo = NULL, xinfo = NULL, weighted = NULL)
#
# r.hub<- rankHub(adj, bdeg = NULL, level = c("x", "y"))
# r.hub[[1]]
# r.hub[[10]]
# r.hub[[15]]Running glasso
load("data_retro.RData")
data<- t(data)
p<- dim(data)[2] # no. of nodes
n=nrow(data)
s<- var(data)
a<-glasso(s, rho=0.05)
a$w[1:10,1:10]## [,1] [,2] [,3] [,4] [,5] [,6]
## [1,] 0.08227358 0.00000000 0.0000000 0.00000000 0.00000000 0.000000e+00
## [2,] 0.00000000 0.09158301 0.0000000 0.00000000 0.00000000 0.000000e+00
## [3,] 0.00000000 0.00000000 0.1060627 0.00000000 0.00000000 0.000000e+00
## [4,] 0.00000000 0.00000000 0.0000000 0.08583836 0.00000000 0.000000e+00
## [5,] 0.00000000 0.00000000 0.0000000 0.00000000 0.08302805 0.000000e+00
## [6,] 0.00000000 0.00000000 0.0000000 0.00000000 0.00000000 1.308527e-01
## [7,] 0.00000000 0.00000000 0.0000000 0.00000000 0.00000000 1.030822e-04
## [8,] 0.00000000 0.00000000 0.0000000 0.00000000 0.00000000 0.000000e+00
## [9,] 0.00000000 0.00000000 0.0000000 0.00000000 0.00000000 3.696276e-05
## [10,] 0.00000000 0.00000000 0.0000000 0.00000000 0.00000000 0.000000e+00
## [,7] [,8] [,9] [,10]
## [1,] 0.0000000000 0.00000000 0.000000e+00 0.00000000
## [2,] 0.0000000000 0.00000000 0.000000e+00 0.00000000
## [3,] 0.0000000000 0.00000000 0.000000e+00 0.00000000
## [4,] 0.0000000000 0.00000000 0.000000e+00 0.00000000
## [5,] 0.0000000000 0.00000000 0.000000e+00 0.00000000
## [6,] 0.0001030822 0.00000000 3.696276e-05 0.00000000
## [7,] 0.1217530179 0.00000000 4.861388e-04 0.00000000
## [8,] 0.0000000000 0.09427377 0.000000e+00 0.00000000
## [9,] 0.0004861388 0.00000000 1.271970e-01 0.00000000
## [10,] 0.0000000000 0.00000000 0.000000e+00 0.09151835
## [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8]
## [1,] 12.15457 0.00000 0.000000 0.0000 0.00000 0.000000 0.000000 0.0000
## [2,] 0.00000 10.91906 0.000000 0.0000 0.00000 0.000000 0.000000 0.0000
## [3,] 0.00000 0.00000 9.428381 0.0000 0.00000 0.000000 0.000000 0.0000
## [4,] 0.00000 0.00000 0.000000 11.6498 0.00000 0.000000 0.000000 0.0000
## [5,] 0.00000 0.00000 0.000000 0.0000 12.04412 0.000000 0.000000 0.0000
## [6,] 0.00000 0.00000 0.000000 0.0000 0.00000 7.642636 0.000000 0.0000
## [7,] 0.00000 0.00000 0.000000 0.0000 0.00000 0.000000 8.305836 0.0000
## [8,] 0.00000 0.00000 0.000000 0.0000 0.00000 0.000000 0.000000 10.6074
## [9,] 0.00000 0.00000 0.000000 0.0000 0.00000 0.000000 0.000000 0.0000
## [10,] 0.00000 0.00000 0.000000 0.0000 0.00000 0.000000 0.000000 0.0000
## [,9] [,10]
## [1,] 0.000000 0.00000
## [2,] 0.000000 0.00000
## [3,] 0.000000 0.00000
## [4,] 0.000000 0.00000
## [5,] 0.000000 0.00000
## [6,] 0.000000 0.00000
## [7,] 0.000000 0.00000
## [8,] 0.000000 0.00000
## [9,] 7.873255 0.00000
## [10,] 0.000000 10.92677
## [1] 205
## degree
glasso.degree=apply(fit.gl, 2, sum)
tab<- cbind.data.frame("Nodes"=colnames(data), "degree"=glasso.degree)
ss<- sort(tab$degree, ind=T, decreasing = T)
tab.sort<- tab[ss$ix,]
tab.sort$Nodes<- factor(tab.sort$Nodes, levels=unique(tab.sort$Nodes))
dat<- tab.sort[1:30,]
#ggplot(dat, aes(y=degree, x=Nodes)) +geom_bar(position="stack", stat="identity") + theme_bw() + theme(axis.text.x = element_text(colour="black",size=15,angle=45,hjust=1,vjust=1,face="plain"),axis.title.x = element_text(colour="black",size=25,angle=0,hjust=.5,vjust=.5,face="plain"), axis.text.y = element_text(colour="black",size=15,angle=0,hjust=1,vjust=0,face="plain"), axis.title.y = element_text(colour="black",size=25,angle=90,hjust=.5,vjust=.5,face="plain"),legend.title = element_text(colour="black",size=15), legend.text=element_text(colour="black",size=15), legend.direction = "horizontal", legend.position = "top") + theme(plot.margin = unit(c(1,1,1,1), "cm"))Running space.neighbor
load("data_retro.RData")
data<- t(data)
p<- dim(data)[2] # no. of nodes
n=nrow(data)
alpha=1
l1=1/sqrt(n)*qnorm(1-alpha/(2*p^2))
iter=10
result1=space.neighbor(data, lam1=l1*0.7, lam2=0)
result1$ParCor[1:10,1:10]## [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
## [1,] 1 0 0 0 0 0.00000000 0 0 0.00000000 0
## [2,] 0 1 0 0 0 0.00000000 0 0 0.00000000 0
## [3,] 0 0 1 0 0 0.00000000 0 0 0.00000000 0
## [4,] 0 0 0 1 0 0.00000000 0 0 0.00000000 0
## [5,] 0 0 0 0 1 0.00000000 0 0 0.00000000 0
## [6,] 0 0 0 0 0 1.00000000 0 0 0.06109301 0
## [7,] 0 0 0 0 0 0.00000000 1 0 0.00000000 0
## [8,] 0 0 0 0 0 0.00000000 0 1 0.00000000 0
## [9,] 0 0 0 0 0 0.06109301 0 0 1.00000000 0
## [10,] 0 0 0 0 0 0.00000000 0 0 0.00000000 1
## [1] 396
## degree
spacenb.degree=apply(fit.space.nb, 2, sum)
tab<- cbind.data.frame("Nodes"=colnames(data), "degree"=spacenb.degree)
ss<- sort(tab$degree, ind=T, decreasing = T)
tab.sort<- tab[ss$ix,]
tab.sort$Nodes<- factor(tab.sort$Nodes, levels=unique(tab.sort$Nodes))
dat<- tab.sort[1:30,]
#ggplot(dat, aes(y=degree, x=Nodes)) +geom_bar(position="stack", stat="identity") + theme_bw() + theme(axis.text.x = element_text(colour="black",size=15,angle=45,hjust=1,vjust=1,face="plain"),axis.title.x = element_text(colour="black",size=25,angle=0,hjust=.5,vjust=.5,face="plain"), axis.text.y = element_text(colour="black",size=15,angle=0,hjust=1,vjust=0,face="plain"), axis.title.y = element_text(colour="black",size=25,angle=90,hjust=.5,vjust=.5,face="plain"),legend.title = element_text(colour="black",size=15), legend.text=element_text(colour="black",size=15), legend.direction = "horizontal", legend.position = "top") + theme(plot.margin = unit(c(1,1,1,1), "cm"))Running space.neighbor
load("data_retro.RData")
data<- t(data)
p<- dim(data)[2] # no. of nodes
n=nrow(data)
alpha=1
l1=1/sqrt(n)*qnorm(1-alpha/(2*p^2))
iter=10
result2=space.joint(data, lam1=6, lam2=0, iter=iter)## [1] "iter=1"
## [1] "iter=2"
## [1] "iter=3"
## [1] "iter=4"
## [1] "iter=5"
## [1] "iter=6"
## [1] "iter=7"
## [1] "iter=8"
## [1] "iter=9"
## [1] "iter=10"
## [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
## [1,] 1 0 0 0 0 0.00000000 0 0 0.00000000 0
## [2,] 0 1 0 0 0 0.00000000 0 0 0.00000000 0
## [3,] 0 0 1 0 0 0.00000000 0 0 0.00000000 0
## [4,] 0 0 0 1 0 0.00000000 0 0 0.00000000 0
## [5,] 0 0 0 0 1 0.00000000 0 0 0.00000000 0
## [6,] 0 0 0 0 0 1.00000000 0 0 0.01794755 0
## [7,] 0 0 0 0 0 0.00000000 1 0 0.00000000 0
## [8,] 0 0 0 0 0 0.00000000 0 1 0.00000000 0
## [9,] 0 0 0 0 0 0.01794755 0 0 1.00000000 0
## [10,] 0 0 0 0 0 0.00000000 0 0 0.00000000 1
## [1] 228
## degree
spacejt.degree=apply(fit.space.jt, 2, sum)
tab<- cbind.data.frame("Nodes"=colnames(data), "degree"=spacejt.degree)
ss<- sort(tab$degree, ind=T, decreasing = T)
tab.sort<- tab[ss$ix,]
tab.sort$Nodes<- factor(tab.sort$Nodes, levels=unique(tab.sort$Nodes))
dat<- tab.sort[1:30,]
#ggplot(dat, aes(y=degree, x=Nodes)) +geom_bar(position="stack", stat="identity") + theme_bw() + theme(axis.text.x = element_text(colour="black",size=15,angle=45,hjust=1,vjust=1,face="plain"),axis.title.x = element_text(colour="black",size=25,angle=0,hjust=.5,vjust=.5,face="plain"), axis.text.y = element_text(colour="black",size=15,angle=0,hjust=1,vjust=0,face="plain"), axis.title.y = element_text(colour="black",size=25,angle=90,hjust=.5,vjust=.5,face="plain"),legend.title = element_text(colour="black",size=15), legend.text=element_text(colour="black",size=15), legend.direction = "horizontal", legend.position = "top") + theme(plot.margin = unit(c(1,1,1,1), "cm"))